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We investigate the dynamical behavior of finite rings of classical spin vectors interacting via 
nearest-neighbor isotropic exchange in an external magnetic field. Our approach is to utilize the 
solutions of a continuum version of the discrete spin equations of motion (EOM) which we derive 
by assuming continuous modulations of spin wave solutions of the EOM for discrete spins. This 
continuous EOM reduces to the Landau-Lifshitz equation in a particular limiting regime. The 
usefulness of the continuum EOM is demonstrated by the fact that the time-evolved numerical 
solutions of the discrete spin EOM closely track the corresponding time-evolved solutions of the 
continuum equation. Of special interest, our continuum EOM possesses soliton solutions, and we 
find that these characteristics are also exhibited by the corresponding solutions of the discrete EOM. 
The robustness of solitons is demonstrated by considering cases where initial states are truncated 
versions of soliton states and by numerical simulations of the discrete EOM equations when the 
spins are coupled to a heat bath at finite temperatures. 
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I. INTRODUCTION 

With the remarkable progress in recent years in syn- 
thesizing and analyzing magnetic molecules of great di- 
versity it is timely to focus on the dynamical be- 
havior of finite arrays of interacting spins. In this paper 
we direct our attention to ring structures consisting of N 
equally-spaced (spacing a) classical spins (unit vectors) 
that interact via ferromagnetic nearest-neighbor isotropic 
exchange and are subject to an external magnetic field 
B = Be. The equations of motion (EOM) of the discrete 
spin vectors can be written as 



— s„(t) = s n (t) x (s n _i(i) 
at 



l(*) 



(1) 



0, 



,N-1, 



where all variables are dimensionless and the unit vectors 
s n are subject to the cyclic condition s n _|_jv — s n . It 
is well known that the solutions of Eq. ([!]) exhibit a 
wide variety of dynamical characteristics. Among these 
are: Exact analytical solutions, called "spin waves" 
numerical solutions that exhibit chaotic behavior as 
well as numerical solutions that are discretized versions 
of the soliton solutions [fj of a continuum version of 
Eq. (U), specifically the Landau-Lifshitz equation [7fl in 
the case of B ^ 0. 



The key result of this paper is that the EOM of 
Eq. |T]) admit a wide class of numerical solutions with 
soliton-like characteristics. They are discretized versions 
of soliton solutions of a new continuum EOM that 
includes the Landau-Lifshitz equation as a special, 
limiting case. We will refer to such discretized solutions 
as "quasi-solitons" , using this cautious terminology since 
we do not know whether they correspond to strict soliton 
solutions of Eq. (TTJ). A single quasi-soliton solution is a 
localized disturbance that propagates in the ring while 
for the most part maintaining its initial shape. We 
also find that it is possible to achieve two-quasi-soliton 
solutions, where a pair of initially separated single 
quasi-solitons collide multiple times and emerge from 
collisions without appreciable modification. Moreover, 
we identify solutions of Eq. (TTJ) that exhibit a transition 
from laminar-like to turbulent-like behavior. 

Our continuum EOM is obtained as follows. We recall 
the general spin wave solution of Eq. (TTJ) 



where 



VT — z 1 cos(qn — Bt — tut) 
s n (t) = | Vl - z 2 sin(qn - Bt - wt) 

z 



9~ t (k = 0,...,N-l) 



(2) 



(3) 



the angular frequency uj is given by 
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U) = 2(1 — COS q)z 



(4) 
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and z is an arbitrary real number in the interval — 1 < 
z < 1. In Eq. $2$ we replace z and u>t by two sets of site- 
dependent functions z n (t) and (p n (t), respectively which 
are assumed to vary slowly over a distance a. As shown 
in Sec. UH in a second order approximation scheme we 
obtain an EOM of the form 



dS 
~dt 



B + (l -cosq) (e- (2S - 

„ . , „, dS 
-2a sing (e • b) 

ax 



S x e 



-a cos q S x 



d 2 S 



(5) 



In the limit q — > 0, on rescaling the coordinate a;, this 
equation reduces to the Landau-Lifshitz equation, 



as 

~dt 



S x 



d 2 S 

dx 2 



B 



(6) 



and hence Eq. can be considered as the generaliza- 
tion of Eq. |6]). As we show below Eq. §5§ possesses 
one-soliton solutions, even if B equals 0, and probably 
also A/"-soliton solutions. 



N = 100. Were we not guided by the solutions of the 
continuum equation, the chance of successfully generat- 
ing a quasi-soliton solution to the discrete EOM is much 
like the proverbial "finding a needle in the haystack" . 
In support of this, we find that for virtually all other 
choices of an initially localized spin configuration the 
spin pattern spreads with time and does not reassemble, 
let alone maintain its localized shape. This also applies 
to the first order approximation of the EOM, where we 
have also found numerical solutions that closely follow 
the initial stages of the evolution towards a caustic. 

It is remarkable that the dynamical behavior of 
a relatively small number of interacting spins is so 
multi-faceted and unexpectedly rich, and that one is 
dependent on a study of the continuum version of the 
EOM to provide the key for unraveling and elucidating 
the dynamical properties of the discrete EOM. We have 
not attempted to analyze the third order version of 
the continuum EOM due to the severe mathematical 
complications that arise. Other interesting dynamical 
behavior, as yet unknown to us, might still be found 
from a deeper mathematical analysis of the complete 
continuum EOM. 



The rationale for our strategy of developing Eq. (JSJ) 
is as follows: The replacement of the finite difference 
expression s, i+ i +s„_i — 2s„ in Eq. |T]) by the single term 
§pr, as in Eq. ©, is accurate as long as the continuous 
function S(x,t) varies sufficiently slowly over a distance 
a so that the effects of an infinite number of higher 
order spatial derivatives can be neglected. In particular, 
the solutions of the Landau-Lifshitz equation can only 
describe small deviations from spin wave solutions with 
q approximately zero. By contrast Eq. ([5]) starts from 
the expressions of Eq. ([2]) for any choice of q allowed by 
Eq. ©. 

A first order version of the continuum EOM of 
Eq. ([5]), as shown in Sec. IIIIl can be solved completely. 
It solutions do not describe solitons but deformations 
of the initial profile similar to those of the Hamilton- 
Jacobi equation for free particles. After a finite time 
caustics are necessarily formed indicating that the 
first order approximation breaks down. In our nu- 
merical studies of Eq. (TT]) we do find instances where 
smooth or "laminar" solutions give way to a chaotic 
or "turbulent" regime, but this does not occur in general. 

Solitons appear in the second order approximation 
of the EOM of Eq. (J5]). Using a discretized version of 
the initial form of a solution of Eq. ([5]), we find that, 
under certain conditions, the time-evolved numerical 
solution of the discrete spin EOM closely tracks the 
time-evolved solution. We call these numerical solutions 
"quasi-solitons" . However, we have found that, in the 
case of B = 0, the existence of quasi-solitons is limited 
to a small "window" of q- values, say ^ < \q\ < for 



Another important aspect is the question of how 
robust quasi-solitons are. One possible perturbation 
of solitons is the coupling of the spins to a heat bath. 
This kind of thermal perturbation is incorporated in our 
numerical treatment of the discrete EOM in Sec. IVII 
Likewise, we explore rings with relatively small values 
of N. We find that quasi-soliton solutions are relatively 
robust, surviving as long-lived spatially localized pat- 
terns as the temperature is increased and/or the initial 
profile is truncated or as N is decreased to as small as 11. 

The kind of quasi-solitons investigated in the present 
paper represent only a certain fraction of solitary 
phenomena in spin systems known from the literature. 
Another class of solitons exists for classical chains 
with a particular form of anisotropy that leads to the 
Sine-Gordon equation 8j. Other soliton solutions for 
continuous EOM describing similar systems have been 
given in Solitary solutions of a quantum spin 

chain have been investigated, for example, in A 
very large literature (see Q for a comprehensive review) 
exists which has shown that these EOM are of direct 
relevance to a wide class of one-dimensional magnetic 
materials and provide a variety of theoretical models in 
quantum field theory 

The layout of this paper is as follows. In Sec. [Til 
starting from the discrete EOM for the interacting 
spin system, Eq. (TT]), we review the well-known exact 
spin-wave solutions, and then proceed to derive a 
continuum version of the EOM describing modulated 
spin waves. Sec. IIIIl is devoted to deriving the full set 
of solutions of the first order version of the continuum 



3 



EOM, and in particular to showing that any solution 
necessarily develops a caustic in the course of time. 
We include examples of numerical solutions of Eq. (fTJ) 
showing such behavior. In Sec. IIVI we give a detailed 
analysis of the soliton solutions of the second order 
version of the continuum EOM. We obtain a formula 
for the dependence of the soliton amplitude A on 
the wave number q, the soliton velocity u, and the 
magnetic field B. In the limit q — > this formula 
reduces to the well-known result for the Landau-Lifshitz 
equation. In Sec. [V] we present results of our numerical 
analysis of the discrete EOM of Eq. ([1]). These include 
one-quasi-solitons and multiple-quasi-solitons starting 
from initial profiles obtained from soliton solutions of 
the second order continuum EOM. Sec. I VII is devoted 
to the effects of starting from initial spin configurations 
that are truncated versions of soliton states, as well as 
including the coupling of the discrete spins to a heat 
bath. Finally, in Sec. IVIll we summarize our results and 
discuss open questions. 



II. MODULATED SPIN WAVES 

We consider classical spin rings with N spins inter- 
acting via ferromagnetic nearest-neighbor coupling. The 
equations of motion are given by Eq. |T]) where B = Be 
is the dimensionless magnetic field and e is chosen as a 
unit vector in the direction of the 3-axis. Mathemati- 
cally, the case of anti-ferromagnetic coupling is included 
since it can be achieved by the transformation s n i— > — s n . 
However, solitons which are excitations of a background 
of fully aligned spins will be thermodynamically stable 
only in the ferromagnetic case. Hence we will stick to 
the ferromagnetic case for physical reasons. For N > 4 
the general solution of Eq. flj) cannot be given in closed 
form. However, as mentioned above, there exist spin 
wave solutions of the form of Eq. @. These solutions are 
parameterized by a variable z, the common 3-component 
of all spin vectors, and by the wave number q which 
assumes the discrete values of Eq. ([3]) . 

In this article we will consider modulations of Eq. 
where z and uit are replaced by two sets of functions z n (t) 
and ip n (t), respectively, which are assumed to vary slowly 
over a distance a. That is, s n (t) is given by 



y /j - z n (t) 2 cos(qn - Bt- <p n (t)) 
s„(i) = ( y/1 - z n {tf sm(qn -Bt- ip n (t)) | (7) 

Zn(t) 

It is convenient to remove the term qn by the transfor- 
mation 



s n (t) ee R 3 (-qn) s n (t) , 



(8) 



where i?a(a) denotes the matrix of a rotation about 
the 3-axis with an angle a. One could also remove the 



magnetic field by a suitable uniform rotation but we 
prefer to retain B in what follows in order to make 
the transition to the Landau-Lifshitz equation more 
transparent. 

In the continuum approximation we introduce smooth 
functions z(x, t), <f(x, t) and S(x,i) such that z(na,t) = 
z n (t) and <p(na,t) — tp n (t), and 



y/l - z(x,t) 2 cos(Bt + ip(x, t)) 
S{x,t) = s x/a (t)= [ -y/l - z(x,t) 2 sm{Bt + tp(x,t)) 
z(x, t) 



(9) 



In the continuum approximation we allow for the pa- 
rameter q to take on any values in the interval (— tt,tv), 
whereas the Landau-Lifshitz equation is obtained in the 
limit q — ► 0. 

Using the addition theorems for cos and sin we obtain 

^ * /l — '/ 2 . - [rnQ rv i rT)Q rt zn. «in n i Qin ri\ \ 



Sfi±l 



where 



[cos a± cos q =F sin a± sin q\ 



1 — %n±\ [ sm a ± cos 1 ^ cos a ± sm l] 

V Zn±l(t) 



a± = qn-Bt- tp n ±i(t) ■ 



(10) 



(11) 



Applying the rotation Rs(— qn) to s„ + i + s„_i hence 
yields 

R3(-qn)(s n+1 +s„_i) = 

(s n+ i +s„_i)cos<7 
+e x (s n+ i - s„-l) sinq 
+(z n+ i - cosg)e. (12) 

Finally, we use the representation of Eq. @ and ap- 
ply the rotation R 3 (-qn) to the EOM of Eq. (D thus 
obtaining the continuum approximation of the form 



(13) 



where 



H(x, t) = [(S(x + a,t) + S(x-a, t))] cos q 
+e x [S(x + a.t) — S(x — a, t)] sing 
+ [B + (z(x + a,t) + z(x - a, *)) (1 - cos q)] e. (14) 

In the following we consider truncated versions of 
Eqs. (|13p and (TT4"]) . These are obtained by first expand- 
ing S(x± a, t) in powers of a. We obtain an infinite order 
partial differential equation (PDE) of the form of Eq. 



4 



with 



H = 2fs + -S" + ^S( 4 > 
2 4! 

-2c • [ a S' + ^S^+. 



cosq 



sm q 



B + 2[z+ yz" + ^-z< 4 ) + . . . ) (1 - cosg) 



(15) 



where the arguments (x, t) have been suppressed, and we 
abbreviate the spatial derivatives of S and z by a prime 
and the time derivatives by a dot in what follows. In 
the nth order approximation we keep terms containing 
powers of a through a". In this article we consider zeroth, 
first, and second order approximations only. The Oth 
order PDE is given by 



S = (B + 2z(l - cosq))S x e. 



(16) 



This equation is solved by the spin wave solution, Eq. 
but it also has more general solutions. 



III. FIRST ORDER EQUATIONS 

In the first order approximation Eqs. (fl3|) and (fl"5|) 
reduce to 

S = [B + 2(1 - cos q) (e • S)] S x e - 2a sin q (e • S) S' . 

(17) 

Using the representation of Eq. ([9]) we can rewrite this 
equation as 

<p = B + 2z(x,t)((l — cosq) — a sing ip') (18) 
z = — 2asinq z z' . (19) 



First we consider Eq. (fT^|) , which is an autonomous equa- 
tion for z(x,t). Its solution is well-known, see e. g. [141 ] . 
Equation (fT9|) essentially describes the velocity field of a 
system of free particles with a given velocity distribution 
for t = 0. To make this more transparent, consider the 
transformation 



v(x, t) = z{x, 
which transforms Eq. (Til)]) into 



t 



2a sin q 



v + vv' = . 



(20) 



(21) 



Hence v is constant along the lines with slope 4| = v 
whence the above kinematical interpretation follows. Let 
£ = x — ut be the cc-coordinate of the intersection of the 
axis t = and the line through the point (x, t) with 
slope u. Hence u = u(£,0). This yields a parametric 
representation of the graph of the general solution x i— > 
v(x,t) of Eq. (|2ip for fixed t, namely 



«(f>0) - 



(22) 
(23) 




200 



FIG. 1: Caustics in a line field with initial slope v(x,0) 

2 + 4 hlU ^ 100 ^ • 



where v(£, 0) is an arbitrary initial condition for the con- 
sidered PDE. An example of the resulting line field is 
given in Fig. [TJ We see that there are domains, where 
the lines in the line field intersect, bounded by curves 
called "caustics". These domains correspond to points 
where the parametric representation of Eqs. (|2"2")l . (J22J) 
defines a multi-valued function. The earliest time when 
caustics appear is given by 



tr 



9£ 



(24) 



where £o is the point with the largest negative derivative 
of v(£,0). Initial profiles v(£,0) with only positive slope 
do not have caustics for t > but cannot be realized on 
a spin ring with periodic boundary conditions. 

The same analysis applies to the original equation 
Eq. (fT^|) , except that the caustic time t c in Eq. (p24"|) has to 
be divided by 2a sin q. Hence the first order PDE predicts 
that after some time t c the profile of the 3-components 
of the spin ring develops an infinite slope. That is for 
later times the 1st order analysis is no longer valid and 
the equations of motion must be studied by using higher 
order approximations. Remarkably, there exist examples 
where the numerical solution of Eq. |T]) is well described 
by its 1st order approximation for times t with < t < t c . 
Close to the time t = t c the smooth "laminar" time evo- 
lution breaks down and a different regime begins to pen- 
etrate the spin ring. We will call this regime "turbulent" 
for sake of simplicity, but we do not know whether it is 
genuinely chaotic or still regular on a shorter length scale. 
From our numerical results the latter cannot be excluded. 
Figure [2] shows the numerical solution z n (t) of Eq. ([T|) for 
anJV= 100 spin ring with q — 167r/100, together with 
the 1st order PDE solution z(x,t) at the time t = t c . 
The two sets of results are in close agreement, although 
at time t — t c the continuum approximation is no longer 
applicable. 

We remark that equation (fT5)) can also be solved in 
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z(x,t) z(x,t) 




FIG. 2: Numerical solutions z n (t) of the discrete EOM of 
Eq. ([1]) for a spin ring with N = 100 (blue dots) and solution 
of the 1st order continuum approximation z(x,t) (solid red 
curve), where t — nt c /5, n = 0,1, ... ,5 and t c is the caustic 
time defined in Eq. (I24|l . For t > t c the solid curve bends over 
and defines a 3-valued function. 



z(x,t) z(x,t) 




20 40 60 80 100 20 40 60 80 100 



FIG. 3: Comparison of a numerical solution for z„ (t) at fixed 
times t = 0, 12, . . . , 60 of the equations of motion for an N = 
100 spin ring (blue dots) with the numerical solution of its 
2nd order continuum approximation Eqs. (|27p . (|28p (solid red 
curve) at the times t = 0, 12, ... , 48. The turbulent regime 
begins to spread at approximately t = 48. By t — 60 the 
numerical solution breaks down and hence cannot be drawn. 



closed form. Its general solution utilizes the solution 
z(x,t) of Eq. (fl"9|) already obtained and the initial an- 
gle distribution ipo(x) — f(x, 0) and it is given by 



<p(x, t) = ipo(x—2a sing z(x, t)t)+z(x, t) 2(1— cos q)t+Bt . 

(25) 

Hence the effects of the caustics will also be manifested in 
the behavior of the azimuthal angles of the spin vectors. 



Substituting Eq. §§§ leads to 



ip = B + 2z (1 — cos q — a sing tp') 

II ( /9 ZZ 

z" + z cos q ip u + ^ 

V 1 — z A 



—2a sinq z z' 

+a 2 cos<j(2z z' tp' -(1- z 2 ) ip"). 



(1 - z 2 ) 2 



(27) 
(28) 



IV. THE SECOND ORDER EQUATIONS 

In the preceding section we have seen that the 1st or- 
der approximation predicts its own break-down since it 
will necessarily lead to a caustic or "turbulent" behavior. 
Hence the question arises whether higher order approx- 
imations give rise to the same or different effects. The 
second order equation is given by 

S = [£ + 2(l-cosg)(e-S)]S x e 
— 2asmq(e ■ S)S' 

+ a 2 S x [cosgS" + (l-cosg)(e-S)e]. (26) 



In solving Eq. (I26p numerically and analytically we 
have found instances where the onset of turbulence is 
accelerated as well as where it is entirely suppressed. We 
do not, however, have a clear understanding about the 
underlying causes for either behavior. As an example 
of the first case, acceleration of turbulence, in Fig. [3] we 
compare results of the numerical solution of the EOM of 
Eq. ([1]) with the numerical solution of Eqs. ([27]) and (j2"5j) . 
The time of the onset of the turbulence is approximately 
2.75 times shorter than that predicted by the 1st order 
approximation. 

The other case, preservation of laminar behavior, 
due to 2nd order effects, can be most impressively 
demonstrated by the existence of solitons, presented in 
the next subsection. 
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A. Exact single soliton solutions 

The mathematical derivation of soliton solutions is to 
a large extent analogous to the Landau-Lifshitz case, see 
@. To simplify the derivation we set a = 1 in the 
following. 

We seek solutions of Eqs. (f2"T|) and (j2"5|) of the form 



(29) 
(30) 



tp(x, t) — <&(x — ut) 
z(x, t) — Z(x — ut) 



The resulting equation for Z{x) can be integrated once 
yielding 



$' = sec q 



u(Z — Zq) — Z 2 smq 

l^Z 2 



(31) 



where Zq is some integration constant. Substituting this 
into the equation resulting from Eq. (f2"T)l yields a sec- 
ond order equation for Z(x). Multiplying with Z'(x) 
and integrating with respect to x yields an expression for 
Z'(x) 2 containing an integration constant c. This expres- 
sion is analogous to a one-dimensional potential if Z' (x) 2 
is viewed as the kinetic energy and the total energy is set 
equal to 0, namely 

Z' 2 = -V(Z) 



2(1- Z 2 (l - cos?)) 
{secg [-1 - 3Z 2 - 2u 2 (2ZZ - 1 - Z$) 

-(1 + Z 2 ) cos2q + 8Z 2 (1 - Z 2 ) sin 4 | ] 

+4 [c(Z 2 - 1) + Z(Z + B(l - Z 2 )) 

+u(Z - Z 3 )t&nq] } . (32) 

V(z) will be called the pseudo-potential. Since — V(Z) is 
biquadratic in Z in the numerator and quadratic in the 
denominator the integral 



x — Xq 



dZ 



(33) 



can be expressed in terms of elliptic functions of the first 
and third kind, see [TtJ • The inverse function of Eq. (|33|) . 
together with the integral of Eq. (|3ip . gives the inverse 
function of the soliton solution. However, the analytic 
expression is very complicated and it is far simpler to 
calculate the integral of Eq. p3|) numerically. The poles 
of V{Z) lie at Z„ = ,, ±:L and hence are outside the 

v y r Vl — cosq 

physical domain Z € [—1, 1] for < |g| < ^. Hence we 
will assume 



< 



7T 

<2 



(34) 



in the following. 



The remaining task is to determine values of the 
parameters q, B,u, Zq,c for which solitary solutions of 
Eqs. (|27p and (|2"5|) exist. The inverse function of 



has the form of a soliton profile if V(Z) has a double 
root at Z = Z\ and a simple root at Z — Z 2 such that 
— 1 < Z\ < Z2 < 1 and V(Z) is negative in the interval 
(Zi, Z2). For sake of simplicity we only consider the case 
Z x = -1. 

The condition that V(Z) has a double root at Z = 
Z\ = — 1 gives us two parameters as a function of the 
remaining two: 



l-B- 
sinq 



(2 + Zq) 2 cosq + (Z 2 - 2) sec<? 
4cosg (l + Z ) 2 



1 



A) 



(35) 
(36) 



The soliton solutions hence depend on three parameters 
which we choose to be B, q and u. A typical form of 
V{Z) is shown in Fig. 2] For the numerical studies in 
Sec. fVl and IVII we mostly chose negative wave-numbers q 
in order to get positive velocities according to Eq. (|36|) . 
The remaining roots Z 2t 3 are 



J 2,3 



where 



- esc 4 — [3 + cos 2q + 2u sin q 
12 

-2cosq(2- BTf(B : q,u))} (37) 



f(B,q,u) = y/B 2 - 2u 2 + 2u sec q(u + B sin q) .(38) 

The condition that V (Z) is negative in the interval Z £ 
(—1, Z 2 ) or, equivalently, V"{—1) < 0, which is necessary 
for the existence of solitons, leads to the inequalities 



2 sin q — 2 cos q\/ (2 — B) sec q — 2 



<u< 



(39) 



2 sin q + 2 cos q\J (2 — B) sec q — 2 



This defines a certain physical domain in the (B, q, u)- 
plane, see Fig. [5l where solitons exist. It is bounded by 
B < 2. The inequalities further imply that only the root 
Z 2 of V(Z) with the —sign in Eq. ((37]) lies inside the 
physical domain (—1, 1). 

Eq. (|3~T|) yields an explicit relation between the ampli- 
tude A = Z 2 — Z\ = 1 + Z 2 and the velocity u of the 
soliton of the form 

A = 2 + using + (B + f(B,q,u))cosq 
4 sin 4 I 

with f(B,q,u) defined by Eq. ([38]) . The two limiting 
forms are 



A = 2+— for<j->0and£<0, (41) 
2B 



to be compared to Eq. (4.7) in Q, and 
.4 



2 + — esc H fusing — |u|-y/2(sccq — l)cosg 
for B -> . 



(42) 
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FIG. 4: The pseudo-potential V{Z) according to Eq. (|32| 
corresponding to the parameters B — 0, q — — n/4 and u = 1 





FIG. 5: The soliton velocities u bounded by a function of q 
and B according to Eq. (|39p . Curves shown correspond to 
B — —2, —1.5, . . . , 1.5 The oo-shaped red curve corresponds 
to B = 0. For B — > 2 the curves shrink to two points at 
q = ±7r/2, w = =F2. 

The physical part of the surface in (g, u, A)-space 
spanned by the lines defined by Eq. ([42]) with q — const, is 
shown in Fig. [5] From these equations it follows that the 
amplitude A diverges for B, q — > and hence no solitons 
exist in this case in accordance with [a] . 



V. NUMERICAL SOLITON SOLUTIONS 

In this section we summarize some of our results for 
numerical solutions of the EOM, Eq. fTJ), in the case 
B = 0, for finite spin rings. The initial spin configu- 
rations arc selected to be soliton profiles as calculated 
from the 2nd order continuum PDE using Eqs. ([31]) 
- (l33l) . Fig. [7] shows the initial profile of a discrete 
soliton for a spin ring with N — 100 and parameters 
q = — 87t/100 and u — — sing. The soliton appears as 
a localized deviation from an otherwise ferromagnetic 
alignment of spins. Using this initial profile we solve 



FIG. 6: Soliton amplitude A as a function of wave number q 
and velocity U for B — according to Eq. (|42() . 



the EOM, Eq. {T]), numerically to determine the time 
evolution of the spin vectors. The results are shown in 
Fig. [51 In panel a) we provide a snapshot of z n (t) for a 
particular time t = 2000. A comprehensive picture for 
all calculated times is given in the color contour plot 
of panel b), where the color gives the value of z n (t) 
according to the coding defined in the legend. The strict 
integrity of the solution as a function of time is striking. 
Throughout the remainder of this paper we use color 
contour plots since they provide a very effective tool for 
visualizing the overall time evolution of the spin profile 
as measured by z n (t). In particular, these plots allow 
one to detect even small changes in the profile, its speed, 
etc. 

Using the initial profile for a single soliton solution it is 
easy to prepare the case of two identical solitons moving 
towards each other. We copy the single soliton profile of 
the N = 100 spin ring onto one half of an N = 200 spin 
ring and its reversed version onto the second half of the 
spin ring. As is seen in Fig. [9l the numerical solution of 
Eq. (dJ describes two solitons that collide and then move 
apart while retaining their initial shapes and velocities. 
This behavior is well-known for exact A/"-soliton solutions 
of other nonlinear wave equations, see, e. g. [l7j . We have 
observed repeated collisions of the soliton pairs without 
noticeable deformation. By using a similar procedure we 
have created the situation where a faster soliton overtakes 
a slower one, and the results are shown in Fig. [10] One 
can see that the penetration of the slower soliton by the 
faster one does not alter the shape or velocity of either 
soliton subsequent to their collision. 
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FIG. 7: Initial profile of a discrete soliton for a spin ring 
with N — 100 calculated from Eqs. (|31f) - (|33[) for parame- 
ters q — — 87r/100 and u = sin |g| . Panel a) shows the in- 
plane view while panel b) shows the structure of the soliton 
as viewed from above. The initial profile z n (t) can be traced 
by eye following the blue arrow tips. The time evolution of 
this soliton is shown in Fig. [8] 



VI. EFFECTS OF TRUNCATION AND FINITE 
TEMPERATURE 

The purpose of this Section is to explore the extent to 
which quasi-solitons are robust to various perturbations. 
We begin by investigating the time evolution of symmet- 
rically truncated versions of initial states that exhibit 
soliton behavior. This is illustrated schematically in 
Fig. [ITJ where the black dots describe a hypothetical 
non-truncated initial profile. In the truncated version, 
for all but M spins the value of z n (0) is set equal to — 1 , 
whereas we retain the original values of z n (0) for the M 
spins. The latter are chosen to be centered about the 
site having the largest value of z n (0). 

Based on the initial profile used for the soliton shown 
in Fig. [5] we have used M = 10 as well as M — 3 for the 
corresponding truncated versions. As shown in Figs. [12] 
and 1131 soliton behavior is clearly maintained for the 
choice M = 10 and surprisingly, despite the fact that 
this is a nonlinear problem, even for M = 3. For the case 
M = 3 there is only an increased level of background 
noise while the major features of the quasi-soliton are 
retained. 

While our quasi-soliton solutions appear to be robust 
against symmetric truncation we have found that 
asymmetric off-center truncation leads to far less stable 
behavior, i. e. , the soliton decays after some time which 
depends on the extent of the asymmetry. 

Another way of truncating the initial profile is to 
reduce the system size N, but leaving the initial profile 
intact. This means, that for the system sizes where 
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FIG. 8: Time evolution of the discrete soliton with initial 
profile shown in Fig. [7] Panel a) shows a snapshot of the 
values of z n (t) for t = 2000. The color coding corresponds to 
the contour plot given in panel b). The data given in panel 
a) correspond to the dashed white line shown in panel b). 



N < M the truncation effects are due to the imposed 
cyclic boundary conditions. Our study is based on the 
numerical solution z n (t) with the initial data calculated 
from Eqs. (f3Tj) - ([33]) using the parameters q = — 327r/100 
and u = sia\q\. In Fig. [H] we show results for a spin 
ring with N — 20. The initial profile is symmetrically 
truncated around the maximum value of z n (t). One 
clearly sees soliton behavior with some background 
noise and a slight pulsation of the amplitude z n {t). 
By reducing the system size to N = 11 spins soliton 
behavior is still visible, however the pulsation of the 
amplitude is strongly enhanced (see Fig. [15|). 

While in all of the above cases the total energy is a 
conserved quantity we now show one specific example of 
the behavior of the quasi-solitons when the spin system 
is coupled to a heat bath. The behavior of classical spin 
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FIG. 9: Collision of two solitons of an N = 200 spin ring. The 
initial profile has been calculated from Eqs. (|3ip - (|33|l using 
the parameters q — ±167r/200 and u — — sin q. The numerical 
solutions are displayed for t = 0, . . . , 4000, including the times 
ti ~ 1000 and ti ~ 2800 where collision between the two 
solitons take place. 
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FIG. 10: Numerical two-soliton solution z n (t) for an N — 200 
spin ring and initial data calculated from Eqs. (13 1 p - (|33[) using 
the parameters q = — 167r/200 and u — 0.625 sin \q\ resp. u = 
1.25 sin \q\. The solution is displayed for t = 0, . . . , 2200. In 
this case two solitons are created which move in the same 
direction while one soliton is twice faster than the other even- 
tually penetrating the slower one. 
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FIG. 11: Non-truncated initial profile (black dots) and trun- 
cated initial profile which is obtained by taking out M spins 
centered about the site having the largest value of z„(0) and 
setting the value of z n (0) equal to —1 for the remaining spins. 
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FIG. 12: Time evolution for an TV = 100 spin ring and an 
initial profile, schematically described by Fig. 1111 where M = 
10 spin vectors are specified as in Fig. [7] and the remaining 90 
spin vectors are initially set anti-parallel to the 3-axis. 



systems in contact with a heat bath can very effectively 
be studied using a constant temperature stochastic spin 
dynamics approach, such as that used in [j| [l6j]. Here 
the spin system is coupled to a heat bath according to a 
Langevin-type approach by including a Landau-Lifshitz 
damping term as well as a fluctuating force with white 
noise characteristics. 

Starting from prescribed initial conditions as those 
used in Fig. [5] the temperature is set to a value T > and 
the trajectory of the system is calculated numerically by 
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FIG. 13: Time evolution for an TV = 100 spin ring and an 
initial profile, schematically described by Fig. 1111 where only 
M — 3 spin vectors are specified as in Fig.[7]and the remaining 
97 spin vectors are initially set parallel to the z-axis. Even 
in this case a rather stable quasi-soliton is created. How- 
ever, with its generation a small localized perturbation ap- 
pears which moves in the opposite direction. Its interference 
with the quasi-soliton can be observed several times (first at 
time t « 100). 



FIG. 15: Time evolution for an TV = 11 spin ring, where 
the initial data are selected as described in the text. While 
some features of soliton-like behavior are still visible, because 
of the small value of TV the results are strongly affected by 
background noise and pulsation of the amplitude z n (t). 
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FIG. 14: Time evolution for an TV = 20 spin ring, where the 
initial data are selected as described in the text. Soliton- 
like behavior persists although some background noise and 
pulsation of the amplitude z n (t) is visible. 
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FIG. 16: Finite temperature effects on a quasi-soliton solution 
Zn(t) for an TV = 100 spin ring and initial data calculated 
from Eqs. (|3ip - (|33[) using the parameters q = — 87r/100 and 
u = sin |g| . Compared to the T = case shown in Fig. [5] 
one observes that thermal fluctuations lead to a decay of the 
soliton, however, it persists for a significant length of time 
before eventually dissipating. 
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solving a stochastic Landau-Lifshitz equation, thereby 
allowing the spin system to exchange energy with its 
environment [j| [HI]. Our results are shown in Fig. [16] 
Thermal fluctuations lead to a decay of the quasi-soliton, 
driving the system towards its (ferromagnetic) equilib- 
rium configuration. However, the soliton persists for a 
significant length of time before eventually dissipating, 
as expected. 

The cases considered in this Section provide a picture 
of quasi-solitons in classical Heisenberg rings that are re- 
markably robust even when perturbed by truncations, 
small system size, and heat bath couplings. 

VII. SUMMARY 

The overriding theme of this paper is that relatively 
small classical Heisenberg rings exhibit a rich variety 
of dynamical properties including soliton behavior. We 
have proposed that an effective strategy for exploring 
the dynamical behavior of these rings is to solve the 
system of discrete EOM of Eq. (TT]) upon selecting initial 
spin configurations that are discretized versions of the 
solutions of a new continuum EOM that we derived 
by allowing for continuous modulations of the spin 
wave solutions given by Eq. In practice, due to 

the complexity of the analysis, we have limited our 
analysis to a second order PDE version, Eq. ([26)) , of the 
continuous EOM of Eqs. (13]) and (HI). This second 
order PDE is a generalization of the Landau-Lifshitz 
equation that has been used extensively in studying 
one-dimensional magnetic systems, e. g. , spin chains. 
Among the class of solutions of Eq. ([25)) we have found 
soliton solutions. Using discretized versions of the 
soliton solutions as input data, we have found that the 
discrete EOM of Eq. (TT]) possess quasi-soliton solutions. 
The robustness of these discrete quasi-solitons has been 
demonstrated both by considering truncated versions of 
continuum soliton solutions as well as by investigating 



the dynamical behavior when the spins are coupled 
to a heat bath at finite temperatures. Apart from 
quasi-soliton behavior there is a vast diversity of other 
dynamical behavior exhibited by finite spin rings. We 
have only sparsely investigated this diversity, and at 
present we have only a limited physical understanding 
of the behavior we have observed. We conjecture that 
the existence of robust quasi-solitons can be explained 
by the existence of exact solitary solutions of Eq. (TJ, 
i. e. , localized solutions satisfying s n+ i(t + T) = s n (t) 
for all n = 1, . . . , N and some fixed T. 

Magnetic solitons have been detected in many mag- 
netic systems through their signature on such observables 
as specific heat, NMR, and neutron scattering Q. Apart 
from any intrinsic interest in the dynamics of classical 
Heisenberg spin rings, we were motivated to undertake 
the present study as a first step in considering their 
quantum analogues. This is especially relevant in view of 
the rapidly growing number of quantum spin rings that 
are now available in the form of synthesized magnetic 
molecules [l[ - @. It is thus timely to ask what is the 
quantum analogue of solitons in ring systems, and how 
might these be detected by experiment [l3|? 
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